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LONG-TERM  GOALS 

The  US  Navy  has  been  developing  its  Multi-Static  Active  Coherent  (MAC)  acoustic  search  system. 
MAC  is  an  air  launched  acoustic  search  system  that  combines  a  newly  developed  coherent  source 
sonobuoy  with  a  field  of  receiver  sonobuoys.  With  the  introduction  of  the  coherent  source,  operational 
search  data  collected  during  MAC  exercises  are  well  suited  for  post-processing  to  infer  waveguide 
properties.  The  long  term  goal  of  this  research  is  to  investigate  the  feasibility  of  estimating  sedient 
acoustic  properties  using  data  collected  during  such  exercises. 

OBJECTIVES 

Modal  dispersion  data  have  been  used  in  inversion  schemes  that  estimate  sediment  properties  in  range- 
independent  and  range-dependent  environments.  The  objective  of  this  work  is  investigate  and 
determine  the  source  and  receiver  characteristics  and  other  experimental  parameters  that  are  required 
to  extract  sediment  properties  with  adequate  resolution  from  mode  dispersion  data  estimated  from  data 
collected  during  operational  activities  of  the  Navy. 

APPROACH 

In  shallow  water,  wave  guide  acoustic  propagation  is  best  described  in  terms  of  the  normal  modes  of 
propagation.  Each  of  these  modes  travel  along  the  waveguide  at  speeds  which  are  dependent  on  the 
acoustic  properties  of  the  waveguide.  The  experimental  set  up  to  extract  sediment  acoustic  properties 
consists  of  source  that  transmits  broadband  signals  and  several  sonobuoy  receiver  systems.  The 
acoustic  field  at  a  receiver  placed  at  a  distance  from  the  source  is  then  collected  and  processed  to 
detennine  the  group  speeds  of  the  modes  at  a  set  of  discrete  frequencies.  This  forms  the  data  set  used 
to  obtain  the  acoustic  properties  of  the  waveguide.  Typically  one  assumes  that  the  water  column 
properties  are  known  by  direct  measurement  and  the  only  unknowns  to  be  estimated  are  the  acoustic 
characteristics  of  the  bottom.  The  approach  for  estimating  the  sediment  acoustic  properties  from 
dispersion  relationship,  i.e.  the  variation  of  the  group  speed  of  each  mode  as  a  function  of  frequency, 
has  been  detailed  in  [1,  2],  This  procedure  is  a  linearized  solution  which  is  based  on  perturbation 
theory.  Consider  a  range-independent  ocean  model  whose  compressional  wave  speed  and  density  are 
represented  by cb(z )  and  pb(z )  ,  respectively.  For  this  model  kn( co)  and  c])„(z,(o)  are  the  eigenvalue 
and  mode  function  of  the  nth  mode,  respectively,  that  satisfy  the  depth-dependent  Helmholtz  equation 
and  boundary  conditions  associated  with  the  waveguide  model.  In  these  expressions  z  refers  to  the 
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depth  and  to  the  frequency  of  the  acoustic  source.  Next,  perturb  the  compressional  wave  speed  by  a 
small  quantity  Ac( z ) .  This  will  result  in  a  change  in  the  group  speed  of  the  propagating  modes.  The 
relationship  between  the  group  speeds  of  the  nth  mode  due  to  a  perturbation  of  the  compressional 
wave  speed  is  given  by  [1] 

— — - - —  =  —  J  — - f  Ac(z)  | (|) t,  (z,  to) | 2  clz.  (1) 

V„ (to)  v, (to)  dm  {  kn (to)  cb  (z) pb (z) 

In  the  above  equation,  vn  and  vn  represent  the  group  speeds  of  mode  n  for  the  unperturbed  and  the 
perturbed  ocean  models,  respectively.  By  modeling  the  sediment  as  a  discrete  set  of  layers  (1)  can  be 
re-written  as  a  matrix  equation  of  the  fonn  Gm  =  d  where  m  represents  the  corrections  to  the 
compressional  wave  speed  values  of  the  sediment  layers,  d  the  data  which  are  the  quantity  on  the  left 
hand  side  of  (1),  and  G  are  obtained  by  converting  the  integral  in  (1)  to  a  sum.  This  matrix  equation  is 
ill  posed  and  special  care  has  to  be  taken  to  solve  it.  A  means  of  regularization  of  the  solution  is 
required  to  obtain  meaningful  solutions  to  the  matrix  equation.  Several  approaches  to  regularization 
such  as  Tikhonov  method  [3]  and  its  extension  Qualitative  Regularization  [4]  have  been  proposed.  It 
has  been  shown  in  [2]  that  the  fonnulation  in  (1)  derived  for  range  independent  environment  can  be 
extended  to  range-dependent  environment  but  it  requires  a  multiplicity  of  source/receiver 
combinations  to  solve  the  problem.  Other  approaches  to  estimating  the  sediment  acoustic  properties 
using  global  optimization  such  as  simulated  annealing  and  genetic  algorithm  have  also  been 
implemented[5,6]. 

WORK  COMPLETED 

The  issues  posed  by  the  use  of  this  approach  for  estimating  the  sediment  acoustic  properties  from  a 
field  of  sonobuoys  together  with  a  compact  source  are: 

•  Method  of  generating  a  signal  of  sufficient  bandwidth 

•  Estimating  the  group  speed  of  modes  using  high  resolution  method 

•  Ability  to  estimate  sediment  properties  when  using  data  from  high  frequencies. 

•  Ability  to  extract  range  dependent  sediment  properties  from  a  field  of  sonobuoys 


A.  Broad  band  signal  generation 

Broad  band  signals  that  have  found  wide  acceptance  are  the  LFM  (Linear  Frequency  Modulated) 
signal  in  which  the  frequency  is  varied  linearly  either  going  up  or  down  with  the  signal  duration. 
Broadband  signals  are  also  created  by  pseudorandom  sequences.  The  maximal  length,  binary,  shift 
register  sequences  are  pseudorandom  sequences  [7],  The  shift  register  sequence  generator  consists  of  a 
shift  register  working  in  conjunction  with  appropriate  logic,  which  feeds  back  a  logical  combination  of 
two  or  more  of  its  stages  to  its  input.  The  output  of  the  sequence  generator  and  the  contents  of  its 
stages  at  any  sample  time  is  a  function  of  the  stages  fed  back  at  the  preceding  sample  time.  If  there  are 
n  stages  in  the  shift  register,  the  maximal  length  sequence  is  a  unique  sequence  of  ones  and  zeros  of 
length  N=2n-\  before  repeating  itself.  For  creating  a  broad  band  signal  the  zeros  in  the  sequence  is 
changed  to  -1,  i.e.  we  now  have  a  sequence  of  1  and  -1  and  then  the  1  are  replaced  by  a  chip  and  the  -1 
are  replaced  by  phase  shifted  version  of  the  chip. 


2 


B.  Simulation  of  data 

Let  us  consider  an  ocean  model  where  the  water  column  is  70  m  in  depth  and  the  sediment  consists  of 
three  layer.  The  thicknesses  of  the  three  layers  are  8  m,  16  m  and  8  m  respectively.  This  is  tenninated 
by  a  half  space  of  constant  compressional  wave  speed  and  density.  The  sound  speed  in  the  water 
column  is  1480  m/s,  the  compressional  wave  speed  in  the  three  layers  are  1569  m/s,  1705  m/s,  and 
1527  m/s  respectively.  The  densities  in  all  the  layers  are  1.6  gm/cc.  The  compressional  wave  speed  in 
the  tenninating  half  space  is  1850  m/s  and  density  1.6  gm/cc.  Attenuation  in  the  sediment  layers  and  in 
the  half  space  have  been  neglected. 

The  data  were  simulated  using  this  ocean  model.  A  broad  band  signal  was  created  using 
pseudorandom  maximum  length  sequence  of  length  51 1. The  sequence  of  1  and  -1  is  used  to  modulate 
a  carrier  signal  with  a  frequency  of  500  Hz  and  each  digit  in  the  sequence  is  chip  of  length  one  cycle  at 
500  Hz.  The  sampling  frequency  was  2500  Hz.  A  time  domain  acoustic  signal  at  the  receiver  can  be 
modeled  by  Fourier  synthesis  of  the  frequency  domain  solution  to  the  acoustic  wave  equation 
p(r,z,co). 

00 

p(r,z;t)  =  J  S(cL>)p(r,z;co)exp(i(L>t)dcL>,  (2) 

-00 

where  S( co)  is  the  spectrum  of  the  transmitted  signal  and  p(r,z;co)  is  the  frequency  domain  solution  at 
range  r  and  depth  z.  To  compute  the  time  domain  solution  above  requires  that  p(r,z;co)  be  calculated 
at  a  number  of  discrete  frequencies,  which  span  the  frequency  band  of  the  source  spectrum.  The  ocean 
environment  was  the  input  to  the  normal  mode  code  KRAKEN  [8]  for  computing  the  frequency 
domain  solution  with  the  range  to  the  receiver,  depths  of  the  source  and  receiver  set  as  required. 

C.  Estimation  of  group  speed  of  modes 

The  group  speeds  as  a  function  of  mode  number  and  frequency  can  be  obtained  from  a  knowledge  of 
the  mode  eigenvalues  and  mode  functions  [9].  This  assumes  that  the  environment  is  fully  known.  In  a 
field  experiment  designed  to  estimate  the  mode  dispersion,  a  broadband  signal  is  broadcast  by  a  source 
placed  in  the  shallow  water  waveguide  and  the  signal  at  a  receiver  placed  at  a  distance  from  the  source 
is  collected  and  processed.  The  group  speed  as  a  function  of  frequency  and  mode  is  obtained  by 
performing  short  time  Fourier  transform  (STFT)  of  the  received  signal.  The  squared  modulus  of  the 
STFT  is  called  the  spectrogram.  Such  a  transform  provides  the  time  for  the  modes  to  travel  from 
source  to  receiver  from  which  the  mode  group  speeds  are  obtained.  If  the  distance  between  the  source 
and  receiver  is  large,  the  modes  are  fully  resolved  and  it  is  not  problematic  to  estimate  the  mode 
dispersion  curves  from  the  spectrogram.  When  the  range  to  receiver  is  short,  the  modes  are  not  well 
resolved  and  a  high  resolution  method  must  be  used  to  resolve  the  modes.  The  high  resolution  method 
proposed  in  [10,  11]  involves  applying  warping  operator  to  the  received  signal. 

In  Figure  1,  two  of  the  five  steps  of  the  warping  procedure  are  shown.  The  top  left  panel  in  the  figure 
is  the  impulse  response  obtained  from  the  received  signal.  The  top  right  panel  shows  the  signal  to 
which  the  warping  operator  has  been  applied.  The  bottom  left  panel  shows  the  spectrogram  of  the 
received  signal.  While  the  lower  order  modes  are  easily  distinguished,  the  higher  order  modes  are  not 
well  resolved.  The  bottom  right  panel  the  shows  time-frequency  representation  of  the  warped  signal. 
The  modes  appear  as  pure  tones  and  are  well  resolved.  The  modes  can  be  filtered  by  simply  masking 
modes  those  that  are  not  of  interest.  The  filtered  modes  are  then  processed  to  estimate  the  mode 
dispersion  data  for  the  filtered  mode.  This  procedure  provides  robust  estimates  of  the  mode  dispersion 
at  frequencies  above  the  Airy  phase.  The  ability  of  inverse  procedure  to  estimate  the  sediment 
properties  in  deeper  layers  depends  on  the  mode  function  having  sufficient  magnitude  at  these  depths. 
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This  requires  inclusion  of  modes  below  the  Airy  phase  frequencies.  Some  improvement  in  this  regard 
is  obtained  by  using  STFT  modified  to  take  into  consideration  the  dispersive  nature  of  the  medium 
[12].  Another  scheme  for  extracting  the  full  dispersion  data  has  been  proposed  [13]  and  this  is  being 
investigated. 
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Fig.  1  The  top  left  panel  is  the  impulse  response  of  the  channel  to  the  transmitted  signal.  The 
topright  panel  shows  the  time  warped  signal.  The  bottom  left  panel  is  the  spectrogram  of  the  impulse 
response.  The  bottom  right  panel  is  the  short  time  Fourier  transform  (STFT)  of  the  warped  signal. 

RESULTS 

The  main  objective  of  this  work  is  to  investigate  the  feasibility  of  extracting  sediment  properties  from 
data  aquired  by  a  distributed  network  of  sonobuoys.  The  source  that  is  likely  to  be  deployed  is  a 
compact  coherent  source.  While  the  details  such  as  the  bandwidth  of  such  source  is  not  known,  the 
tests  for  geoacoustic  inversion  was  done  assuming  that  the  source  will  have  a  bandwidth  of  about  200 
Hz. The  inversions  were  carried  out  with  mode  dispersion  data  in  the  frequency  intervals  50  Hz  -  150 
Hz,  100  Hz  -  300  Hz  and  300  Hz  -  500  Hz. 

A.  Inversion  based  on  perturbation  theory. 

The  inversion  for  the  sediment  compressional  wave  speed  profiles  were  perfored  using  six  different 
sets  of  data.  These  are  (a)  noise  free  data  in  the  frequency  inertvals  50  Hz  -  150  Hz,  100  Hz  -  300  Hz, 
and  300  Hz  -  500  Hz,  (b)  data  with  noise  in  the  frequency  inerval  100  Hz  -  300  Hz,  and  (c)  data  from 
short  range  (5  km  range)  in  the  frequency  inervals  50  Hz  -150  Hz  and  100  Hz  -  300  Hz. 

The  unknowns  to  be  estimated  are  the  compressional  wave  speed  values  in  the  three  layers.  All  the 
other  parameters  were  assumed  known.  The  data  used  in  the  inversion  are  shown  in  Fig  2.  The  results 
of  inversion  are  shown  in  Fig  3.  The  plots  3(a),  3(b),  and  3(c)  in  the  figure  are  the  inversion  results  in 
which  data  were  obtained  from  noise  free  signal  at  the  receivers.  It  is  seen  that  the  inverse  resulted  in 
estimating  the  compressionalwave  speed  in  the  three  layers  without  any  significant  error  when  data 
from  frequency  intervals  50  Hz  -150  Hz  and  100  Hz  -  300  Hz  were  used.  However,  when  data  from 
300  Hz  -500  Hz  were  used  the  inversion  could  not  estimate  the  compressional  wave  speed  in  the  third 
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layer.  In  this  case  the  magnitude  of  the  mode  function  at  depths  corresponding  to  the  deeper  layers  is 
near  zero  and  hence  the  inversion  procedure  was  not  able  to  estimate  the  layer  speed  in  the  third 
sediment  layer. 
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Fig.  2.  The  top  panel  shows  the  estimates  of  group  speed  of  modes  1  to  12  for  frequencies  from  50 
Hz  to  500  Hz.  The  blue  dots  are  the  values  obtained  by  analysis  of  the  noise  free  received  signal 
using  the  warping  method  and  the  red  circle  denotes  the  values  obtained  from  theory.  The  bottom 
left  panel  shows  the  data  for  modes  1  to  7 for  frequencies  from  50  Hz  to  300  Hz  in  the  case  where 
noise  has  been  added  to  the  data.  The  bottom  right  panel  shows  the  data  when  a  range  to  receiver  is 

short  (5  km). 

We  now  examine  the  effect  of  noise  in  the  data.  Gaussian  random  noise  with  SNR  of  5dB  was  added 
to  the  noise  free  signal  at  the  receiver.  The  group  velocities  estimated  from  the  noisy  signal  are 
compared  with  the  true  group  speed  in  the  bottom  left  panel  of  Fig  2.  It  is  seen  that  there  are  a  number 
of  data  points  where  the  estimate  of  group  speeds  from  the  noisy  signal  are  in  error.  In  order  to 
perform  the  inversion  with  this  data,  these  data  points  were  assumed  to  have  larger  variance  in  their 
values  and  this  was  incorporated  in  the  inversion  scheme.  The  result  of  the  inversion  using  the  noisy 
data  is  shown  in  Fig  3(d).  It  is  seen  that  noise  in  the  data  results  in  large  error  in  the  estimation  of  the 
layer  speeds.  The  data  for  the  analysis  discussed  so  far  had  been  obtained  with  the  source  to  receiver 
range  set  at  10  km  or  more.  Since  in  field  experiment  such  a  large  aperture  may  not  be  possible,  data 
was  simulated  with  the  source  to  receiver  range  set  at  5  km.  The  group  speed  estimated  from  the 
simulated  data  with  range  set  at  5  km  and  the  true  speeds  are  shown  in  the  bottom  right  panel  of  Fig  2. 
The  results  of  inversion  using  this  data  are  shown  in  Fig  3(e)  and  3(f).  In  this  case  the  inversions  in  the 
frequency  interval  100  Hz  -  300  Hz  shows  larger  errors  in  the  group  speed  of  the  layers  as  compared 
to  3(b)  which  is  the  case  for  larger  range. 
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Fig.  3.  The  figure  shows  the  staring  model  (blue  firm  line),  the  true  bottom  model  (black  dashed 
line),  and  the  inversion  result  (red  dashed  line).  The  figures  (a),  (b)  and  (c)  use  data  estimated  from 
the  time  frequency  analysis  of  the  noise  free  signal  acquired  at  the  receiver  for  frequencies  50  Hz  - 
150  Hz,  100  Hz-  300  Hz,  and  300  Hz  -500  Hz  respectively.  The  figure  (d)  has  the  result  for  data 
from  noisy  signal  for  frequencies  100  Hz-  300  Hz.  The  figures  (e)  and  (f)  uses  data  from  signal 
acquired  at  a  range  of  5  km  for  frequencies  50  Hz  -  150Hz  and  100  Hz-  300  Hz  respectively. 

B.  Inversion  based  on  global  optimiation. 

Inversion  can  also  be  treated  as  a  parameter  estimation  problem  where  the  attempt  is  to  find  the  bottom 
model  that  predicts  the  measured  quantity  with  minimum  error.  Specifically  let  dk  be  a  vector 
containing  the  data,  i.e.  a  vector  of  the  group  speed  for  a  set  of  frequencies  in  the  case  of  the  kth  mode 
as  detennined  from  field  measurements.  For  a  given  bottom  model  m  which  contains  the  parameters  to 
be  estimated,  the  calculated  group  speed  of  mode  k  is  dk(m).  The  vector  m  contains  the  sediment 
properties  that  are  to  be  estimated.  The  estimation  problem  then  reduces  to  the  detennination  of  the 

parameters  of  m  that  minimizes  Y,i\dk  —  dk(m )  |  .  The  most  commonly  used  global  optimization 
methods  that  have  been  used  to  solve  problems  of  this  kind  are  (a)  Simulated  annealing  and  (b) 

Genetic  algorithm  and  their  variants. 

We  use  the  Genetic  algorithm  to  estimate  the  compressional  wave  speed  in  the  sediment  layers  and  in 
the  half  space.  More  detailed  explanation  of  the  procedure  used  in  the  execution  of  Genetic  Algorithm 
can  be  found  in  [14].  The  search  bounds  and  the  estimates  of  the  parameter  values  are  in  Table  1. 

These  values  were  obtained  using  the  simulated  data  without  and  with  noise  in  the  frequency  interval 
50  Hz  to  150  Hz. 
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Table  1:  The  results  of  inversion  usiiu 

•  genetic  algorithm 

Type  of  data 

Parameter 

Parameter 
Value  (m/s) 

Search  bound 
(m/s) 

Inversion  result 

50  Hz  -  150  Hz 
(No  noise) 

Layer  1 

1569.0 

1500  -  1800 

1575.0 

Layer  2 

1705.0 

1500-  1800 

1714.0 

Layer  3 

1527.0 

1500  -  1800 

1521.0 

Half  space 

1850.0 

1800  -2000 

1852.5 

50  Hz -150  Hz 
(Noise  added) 

Layer  1 

1569.0 

1500  -  1800 

1584.0 

Layer  2 

1705.0 

1500  -  1800 

1710.0 

Layer  3 

1527.0 

1500  -  1800 

1523.0 

Half  space 

1850.0 

1800-  2000 

1918.0 

CURRENT  WORK 

One  of  the  drawbacks  of  the  warping  method  for  extraction  of  the  mode  dispersion  data  is  its  inability 
to  detennine  the  mode  dispersion  at  frequencies  below  that  of  the  Airy  phase  for  that  mode.  A  scheme 
has  been  proposed  [13]  that  will  enable  extraction  of  the  complete  mode  dispersion  data.  The  scheme 
uses  a  tiling  in  the  time-frequency  plane  which  accounts  for  dispersion.  The  change  in  the  tiling 
requires  a  knowledge  of  the  group  velocity  dispersion.  Investigations  on  the  applicability  of  this 
scheme  by  the  use  of  approximate  values  of  the  group  velocity  disperson  is  under  investigation. 

Further  another  approach  suggested  in  [15]  is  also  being  persued.  A  comprehensive  investigation  is 
also  needed  to  study  the  impact  of  imperfect  knowledge  of  the  water  column  properties,  bathymetry 
between  source  and  receiver  and  other  environmental  parameters  on  inversion  for  sediment  properties, 
estimation. 

IMPACT/APPLICATIONS 

The  data  collected  during  this  experiment  will  enable  validation  of  the  proposed  method  for  estimating 
range-dependent  sediment  compressional  wave  speed  from  modal  dispersion  data.  Using  a  distributed 
set  of  receivers  and  a  broadband  source  it  will  be  possible  to  estimate  the  compressional  wave  speed 
profiles  over  a  wide  area.  This  will  therefore  be  a  useful  tool  for  estimating  the  sediment  acoustic 
properties  from  data  collected  during  routine  naval  operations. 

RELATED  PROJECTS 

The  estimation  of  sediment  properties  in  a  shallow  water  environment  is  an  active  area  of  research. 

The  new  procedures  that  are  likely  to  come  out  of  this  work  will  therefore  be  of  great  benefit. 
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